Enteropathogen seroepidemiology among children in low-resource settings

Supplementary Information File 5. Sensitivity analysis of change in IgG used to identify presumed unexposed measurements in Haiti and Kenya.

Notebook summary

Identifying presumed seroconverters

One approach to identifying a seropositivity cutoff in a longitudinal study is to identify children who presumably seroconvert from negative to positive based on a change in their antibody levels, and then use the distribution of their antibody levels before seroconversion to identify a seropositivity cutoff. The rationale for this two-step process that estimates a cutoff (and does not simply classify children as seroconverers versus not) is that the cutoff can then be applied to identify children who were seropositive but did not actually seroconvert during the study. This was common for antibody response to enteric pathogens in the study cohorts, where exposure ocurred frequently and very early in life.

Level of increase in antibody levels to identify seroconversion

The method relies on choosing the appropriate magnitude of change in antibody levels used identify presumed seroconverters, which is generally unknown. Most immunogenicity abnd pathogen challenge studies have used a 4-fold increase. In this study, we used an increase of +2 or more on the log\(_{10}\) scale, equal to a 100-fold increase, to identify presumed seroconverters. Although this was likely a conservative assumption, it still led to high levels of agreement with other classification approachs (generally >95%) when applied across cohorts and antibodies.

Sensitivity analysis description

Here, we conduct a sensitivity analysis that searches over a range of increases in antibody levels used to identify “presumed seroconverters”. For each increase, we then calculated the mean and SD of the distribution before the change, and used that as a seropositivity cutoff. We then compared the classification accuracy of each value against classifications based on cutoffs derived from ROC-curves or mixture models, as available.

In the present study, we used two additional methods to determine seropositivity cutoffs:

• For Giardia (VSP-3, VSP-5), Cryptosporidium (Cp17, Cp23), and E. histolytica (LecA) antigens, known negative and positive samples were used to identify an optimal cutoff based on a receiver-operator characteristic (ROC) curve.
  
• For all other pathogens, we fit a 2-component, finite Gaussian mixture model to the antibody distributions, and then used the lower component’s mean plus 3 standard deviations to define a cutoff value. This identified reasonable cutoff values for many, but not all, additional antibodies. In some cases the mixture model based approach did not identify reasonable cutoff values (i.e., near or above the upper limit assay’s dynamic range).

Given this other classification information, we compared a range of changes in antibody levels used to identify children who presumably seroconverted in terms of agreement with ROC-based and mixture model-based classification of seropositivity.

We searched over a range of changes from +0.3 to +2.5 on the log10 MFI-bg scale (this corresponds to a 10^0.3 = 2-fold to 10^2.5 = 316-fold increase in antibody levels). For each change in antibody levels (0.3 to 2.5), we identified the children who would be identified as seroconverters. From these children, we then estimated the mean plus 3 standard deviations of their antibody levels at their first measurement, when they were presumably unexposed. This is the seropositivity cutoff among the “presumed unexposed”. We then estimated agreement and Cohen’s Kappa statistics for each change and identify the magnitude of change that corresponded with the maximum level of agreement.

Sensitivity analysis results

The optimal change in MFI levels used to identify presumed seroconverters varied by cohort and antigen, where “optimal” was defined as the minimum change in MFI with maximum agreement with respect to ROC-based or mixture model-based seropositivity classifications.

For antigens with both ROC-based cutoffs and mixture model-based cutoffs, results were similar between the two comparisons. This result reflects good agreement between ROC- and mixture model-based cutoffs (detailed results in this article’s Supplementary Information File 2).

Over the range of values considered, the method led to cutoff values with classification agreement that was very high compared with both ROC- and mixture model-based classifications. In Haiti, optimal agreement exceeded 0.92 for all comparisons. In Kenya, optimal agreement exceeded 0.98 for all comparisons.

In Haiti, the change in MFI values that led to highest agreement were generally \(\geq 2\) on the log\(_{10}\) scale. In Kenya, change values were lower and generally in the range of 0.5 to 1.5. In general, across both countries, there was not a large decline in agreement at higher change values imposed to identify presumed seroconverters when developing a cutoff (an exception was Cryptosporidium Cp23 in Kenya).

Taken together, this sensitivity analysis supports the use of an increase of +2 log\(_{10}\) MFI to identify presumed seroconverters in the definition of seropositivity cutoffs for antigens considered.

Script preamble

#-----------------------------
# preamble
#-----------------------------
library(here)

here() starts at /Users/benarnold/enterics-seroepi

here()

[1] "/Users/benarnold/enterics-seroepi"

library(tidyverse)

── Attaching packages ──────────────────────────────────────────────── tidyverse 1.2.1 ──
✔ ggplot2 3.0.0     ✔ purrr   0.2.4
✔ tibble  1.4.2     ✔ dplyr   0.7.4
✔ tidyr   0.8.0     ✔ stringr 1.3.1
✔ readr   1.1.1     ✔ forcats 0.3.0
── Conflicts ─────────────────────────────────────────────────── tidyverse_conflicts() ──
✖ dplyr::filter() masks stats::filter()
✖ dplyr::lag()    masks stats::lag()

library(psych)

Attaching package: ‘psych’

The following objects are masked from ‘package:ggplot2’:

    %+%, alpha

library(gridExtra)

Attaching package: ‘gridExtra’

The following object is masked from ‘package:dplyr’:

    combine

# set up for parallel computing
# configure for a laptop (use only 3 cores)
library(foreach)

Attaching package: ‘foreach’

The following objects are masked from ‘package:purrr’:

    accumulate, when

library(doParallel)

Loading required package: iterators
Loading required package: parallel

registerDoParallel(cores=3)
# bright color blind palette:  https://personal.sron.nl/~pault/ 
cblack <- "#000004FF"
cblue <- "#3366AA"
cteal <- "#11AA99"
cgreen <- "#66AA55"
cchartr <- "#CCCC55"
cmagent <- "#992288"
cred <- "#EE3333"
corange <- "#EEA722"
cyellow <- "#FFEE33"
cgrey <- "#777777"

Leogane, Haiti

Load the formatted data

#-----------------------------
# load the formatted data
# created with 
# haiti-enteric-ab-data-format.Rmd ->
# haiti-enteric-ab-distributions.Rmd
#-----------------------------
dl <- readRDS(here("data","haiti_analysis2.rds"))
# list the enteric antigens and formatted labels for them
mbavars <- c("vsp3","vsp5","cp17","cp23","leca","salb","sald","etec","norogi","norogii")
mbalabs <- c("Giardia VSP-3","Giardia VSP-5","Cryptosporidium Cp17","Cryptosporidium Cp23","E. histolytica LecA","Salmonella LPS group B","Salmonella LPS group D","ETEC LT β subunit","Norovirus GI.4", "Norovirus GII.4.NO")

Sensitivity analysis over fold-change values, ROC-based reference

Giardia, Cryptosporidium, and E. histolytica antigens have ROC-based cutoff values in this study. For these antigens, search over fold-change values to determine the cutoff with best agreement with ROC-based classification.

droc <- dl %>%
  filter(!is.na(roccut)) %>%
  droplevels()
gridroc <- foreach(ab=levels(droc$antigenf),.combine=rbind) %:%
  foreach(icut=seq(0.3,2.5,by=0.01),.combine=rbind) %dopar% {
    
    # downsample to measurements among children <=1 y 
    # who increased by the increment icut in the next period    
    # these are "presumed unexposed" under this cutoff in change icut
    # estimate mean and SD among the presumed unexposed to use as a 
    # seropositivity cutoff
    dmu <- droc %>% 
      ungroup() %>%
      select(antigenf,age,logmfi,logmfidiff) %>%
      filter(antigenf==ab & age <= 1 & logmfidiff>icut) %>%
      summarize(mu=mean(logmfi),sd=sd(logmfi),ucut=mu+3*sd)
    
    # classify all observations using the derived cutoff
    di <- droc %>%  
      filter(antigenf==ab) %>%
      mutate(posi=factor(ifelse(logmfi>dmu$ucut,1,0),levels=c(0,1)))
    
    # estimate agreement
    tabi <- as.vector(table(as.factor(di$posroc),as.factor(di$posi)))
    names(tabi) <- c("roc0cuti0","roc1cuti0","roc0cuti1","roc1cuti1")
    kappai <- psych::cohen.kappa(table(di$posroc,di$posi))$kappa
    tabi <- data.frame(t(tabi))
    tabi <- tabi %>%
      mutate(Nobs=roc0cuti0+roc1cuti0+roc0cuti1+roc1cuti1,
             agreement=(roc0cuti0+roc1cuti1)/Nobs,
             kappa=kappai,
             antigenf=ab,
             abchange=icut,
             abcutoff=dmu$ucut) 
    
  }
gridroc <- gridroc %>% mutate(antigenf=factor(antigenf,levels=levels(droc$antigenf)))
# identify the optimal change based on % agreement, take minimum change
optchgroc <- gridroc %>%
  group_by(antigenf) %>%
  mutate(maxagr=max(agreement,na.rm=TRUE),
         maxchg=ifelse( abs(agreement-maxagr)<0.001,1,0)) %>%
  filter(maxchg==1) %>%
  arrange(antigenf,abchange) %>%
  filter(row_number()==1) %>%
  select(antigenf,optchg=abchange,optagr=agreement,optkap=kappa,optcut=abcutoff)
# merge in the optimal change values
gridroc <- left_join(gridroc,optchgroc,by="antigenf")  %>%
  group_by(antigenf)

Figure to summarize results

General figure scheme:

# this general figure scheme is specific to results stored in
# gridroc and gridmix objects.
plot_optchange <- function(results,classlab="XXX reference standard") {
  
  nantigens <- length(levels(results$antigenf))
  
  # Figure of cutoff by antibody change threshold
  plot_cutoff <- ggplot(data=results,aes(x=abchange)) +
  geom_line(aes(y=abcutoff),color=cmagent) +
  geom_vline(aes(xintercept=optchg),lty="dashed")+ 
  geom_hline(aes(yintercept=optcut),lty="dashed")+ 
  # labels for optimal change and corresponding cutoff based on % agreement
  geom_label(data=filter(results,row_number()==1),
            aes(x=optchg,y=0.5,label=sprintf("%1.2f",optchg)),
            position=position_nudge(x=0),
            fill="white",label.size=0) + 
  geom_label(data=filter(results,row_number()==1),
            aes(y=optcut,x=0,label=sprintf("%1.2f",optcut)),
            position=position_nudge(y=0,x=0.4),
            fill="white",label.size=0) + 
  facet_wrap(~antigenf,nrow=1,ncol=nantigens) +
  # aesthetics
  scale_y_continuous(breaks=seq(0.5,4,by=0.5))+
  scale_x_continuous(breaks=seq(0,3,by=0.5))+
  coord_cartesian(ylim=c(0.5,4),xlim=c(0,3)) +
  labs(x="log10 change in MFI used to identify presumed seroconverters",y="Estimated seropositivity cutoff (log10 MFI-bg)") +
  theme_minimal()
  
  # Figure of Kappa by antibody change threshold
  plot_kappa <- ggplot(data=results,aes(x=abchange)) +
  # geom_line(aes(y=agreement),color=cgreen) +
  geom_line(aes(y=kappa),color=cblue)+
  geom_vline(aes(xintercept=optchg),lty="dashed")+ 
  geom_hline(aes(yintercept=optkap),lty="dashed")+ 
  # labels for optimal change and corresponding kappa based on % agreement
  geom_label(data=filter(results,row_number()==1),
            aes(x=optchg,y=0,label=sprintf("%1.2f",optchg)),
            position=position_nudge(x=0),
            fill="white",label.size=0) + 
  geom_label(data=filter(results,row_number()==1),
            aes(y=optkap,x=0,label=sprintf("%1.2f",optkap)),
            position=position_nudge(y=0,x=0.4),
            fill="white",label.size=0) + 
  facet_wrap(~antigenf,nrow=1,ncol=nantigens) +
  # aesthetics
  scale_y_continuous(breaks=seq(0,1,by=0.1))+
  scale_x_continuous(breaks=seq(0,3,by=0.5))+
  coord_cartesian(ylim=c(0,1),xlim=c(0,3)) +  
  labs(x="",y=paste("Cohen's Kappa,",classlab)) +
  theme_minimal()
  
  # Figure of Agreement by antibody change threshold
  plot_agree <- ggplot(data=results,aes(x=abchange)) +
  geom_line(aes(y=agreement),color=cgreen) +
  # geom_line(aes(y=kappa),color=corange)+
  geom_vline(aes(xintercept=optchg),lty="dashed")+ 
  geom_hline(aes(yintercept=optagr),lty="dashed")+ 
  # labels for optimal change and agreement based on % agreement
  geom_label(data=filter(results,row_number()==1),
            aes(x=optchg,y=0,label=sprintf("%1.2f",optchg)),
            position=position_nudge(x=0),
            fill="white",label.size=0) + 
  geom_label(data=filter(results,row_number()==1),
            aes(y=optagr,x=0,label=sprintf("%1.2f",optagr)),
            position=position_nudge(y=0,x=0.4),
            fill="white",label.size=0) + 
  facet_wrap(~antigenf,nrow=1,ncol=nantigens) +
  # aesthetics
  scale_y_continuous(breaks=seq(0,1,by=0.1))+
  scale_x_continuous(breaks=seq(0,3,by=0.5))+
  coord_cartesian(ylim=c(0,1),xlim=c(0,3)) +
  labs(x="",y=paste("Agreement,",classlab)) +
  theme_minimal()
  
  # stack the figures into a single panel
  grid.arrange(plot_agree,plot_kappa,plot_cutoff)
}

The bottom panel plots the estimated seropositivity cutoff as a function of change in MFI levels used to identify presumed seroconverters. The mean plus 3 SDs of the distribution of antibody levels before presumed seroconversion was used to estimate the cutoff. For a given cutoff, the middle and top panels summarize classification agreement with the ROC-based seropositiity cutoff based on Cohen’s Kappa (middle) and proportion agreement (top).

Labeled dashed lines indicate the minimum magnitude of change in antibody levels that led to highest agreement, along with the corresponding seropositivity cutoff values, Cohen’s Kappa, and agreement.

plot_optchange(results=gridroc,classlab="ROC-based classification")

Sensitivity analysis over fold-change values, mixture model-based reference

Search over fold-change values to determine the cutoff with best agreement with mixture model-based classification

dmix <- dl %>%
  filter(!is.na(mixcut)) %>%
  droplevels()
gridmix <- foreach(ab=levels(dmix$antigenf),.combine=rbind) %:%
  foreach(icut=seq(0.3,2.5,by=0.01),.combine=rbind) %dopar% {
    # downsample to measurements among children <=1 y 
    # who increased by the increment icut in the next period
    # these are "presumed unexposed" under this cutoff in change icut
    # estimate mean and SD among the presumed unexposed to use as a 
    # seropositivity cutoff
    dmu <- dmix %>% 
      ungroup() %>%
      select(antigenf,age,logmfi,logmfidiff) %>%
      filter(antigenf==ab & age <=1 & logmfidiff>icut)
    
    # if there are <2 obs, cannot estimate mean + SD, so skip
    if(nrow(dmu)>=2) {
      
      dmu <- dmu %>%
        summarize(mu=mean(logmfi),sd=sd(logmfi),ucut=mu+3*sd)
    
      # classify all observations using the derived cutoff
      di <- dmix %>%  
      filter(antigenf==ab) %>%
      mutate(posi=factor(ifelse(logmfi>dmu$ucut,1,0),levels=c(0,1)),
             posmix=factor(posmix,levels=c(0,1)))
    
      tabi <- as.vector(table(di$posmix,di$posi))
      names(tabi) <- c("mix0cuti0","mix1cuti0","mix0cuti1","mix1cuti1")
      kappai <- psych::cohen.kappa(table(di$posmix,di$posi))$kappa
      tabi <- data.frame(t(tabi))
      tabi <- tabi %>%
        mutate(Nobs=mix0cuti0+mix1cuti0+mix0cuti1+mix1cuti1,
             agreement=(mix0cuti0+mix1cuti1)/Nobs,
             kappa=kappai,
             antigenf=ab,
             abchange=icut,
             abcutoff=dmu$ucut) 
      
    } else{
      tabi <- data.frame(mix0cuti0=NA,mix1cuti0=NA,mix0cuti1=NA,mix1cuti1=NA,Nobs=NA,agreement=NA,kappa=NA,antigenf=ab,abchange=icut,abcutoff=NA)
      
      }
    }  
gridmix <- gridmix %>% mutate(antigenf=factor(antigenf,levels=levels(dmix$antigenf)))
# identify the optimal change based on % agreement, take minimum change
optchgmix <- gridmix %>%
  group_by(antigenf) %>%
  mutate(maxagr=max(agreement,na.rm=TRUE),
         maxchg=ifelse( abs(agreement-maxagr)<0.001,1,0)) %>%
  filter(maxchg==1) %>%
  # mutate(agr99=ifelse(agreement>0.99 & !is.na(agreement),1,0)) %>%
  # filter(agr99==1) %>%
  arrange(antigenf,abchange) %>%
  filter(row_number()==1) %>%
  select(antigenf,optchg=abchange,optagr=agreement,optkap=kappa,optcut=abcutoff)
gridmix <- left_join(gridmix,optchgmix,by="antigenf") %>%
  group_by(antigenf)

Figure to summarize results

The bottom panel plots the estimated seropositivity cutoff as a function of change in MFI levels used to identify presumed seroconverters. The mean plus 3 SDs of the distribution of antibody levels before presumed seroconversion was used to estimate the cutoff. For a given cutoff, the middle and top panels summarize classification agreement with the mixture model-based seropositiity cutoff based on Cohen’s Kappa (middle) and proportion agreement (top).

Labeled dashed lines indicate the minimum magnitude of change in antibody levels that led to highest agreement, along with the corresponding seropositivity cutoff values, Cohen’s Kappa, and agreement.

plot_optchange(results=gridmix,classlab="Mixture model-based classification")

Asembo, Kenya

Load the formatted data

#-----------------------------
# load the formatted data
# created with 
# asembo-enteric-ab-data-format.Rmd ->
# asembo-enteric-ab-distributions.Rmd
#-----------------------------
dl <- readRDS(here("data","asembo_analysis2.rds"))
# list the enteric antigens and make formatted labels for them
mbavars <- c("vsp3","vsp5","cp17","cp23","leca","salb","sald","etec","cholera","p18","p39")
mbalabs <- c("Giardia VSP-3","Giardia VSP-5","Cryptosporidium Cp17","Cryptosporidium Cp23","E. histolytica LecA","Salmonella LPS group B","Salmonella LPS group D","ETEC LT β subunit","Cholera toxin β subunit","Campylobacter p18","Campylobacter p39")
dl <- dl %>% mutate(antigenf=factor(antigenf,levels=mbalabs))

Sensitivity analysis over fold-change values, ROC-based reference

Giardia and Cryptosporidium antigens have ROC-based cutoff values in this study. For these antigens, search over fold-change values to determine the cutoff with best agreement with ROC-based classification.

droc <- dl %>%
  filter(!is.na(roccut)) %>%
  mutate(antigenf=droplevels(antigenf))
gridroc <- foreach(ab=levels(droc$antigenf),.combine=rbind) %:%
  foreach(icut=seq(0.3,2.5,by=0.01),.combine=rbind) %dopar% {
    
    # downsample to first measurements who increased by the increment icut
    # these are "presumed unexposed" under this cutoff in change icut
    # estimate mean and SD among the presumed unexposed to use as a 
    # seropositivity cutoff
    dmu <- droc %>% 
      ungroup() %>%
      select(antigenf,time,logmfi,logmfidiff) %>%
      filter(antigenf==ab & time=="A" & logmfidiff>icut) %>%
      summarize(mu=mean(logmfi),sd=sd(logmfi),ucut=mu+3*sd)
    
    # classify all observations using the derived cutoff
    di <- droc %>%  
      filter(antigenf==ab) %>%
      mutate(posi=factor(ifelse(logmfi>dmu$ucut,1,0),levels=c(0,1)))
    
    # estimate agreement
    tabi <- as.vector(table(as.factor(di$posroc),as.factor(di$posi)))
    names(tabi) <- c("roc0cuti0","roc1cuti0","roc0cuti1","roc1cuti1")
    kappai <- psych::cohen.kappa(table(di$posroc,di$posi))$kappa
    tabi <- data.frame(t(tabi))
    tabi <- tabi %>%
      mutate(Nobs=roc0cuti0+roc1cuti0+roc0cuti1+roc1cuti1,
             agreement=(roc0cuti0+roc1cuti1)/Nobs,
             kappa=kappai,
             antigenf=ab,
             abchange=icut,
             abcutoff=dmu$ucut) 
    
  }
gridroc <- gridroc %>% mutate(antigenf=factor(antigenf,levels=levels(droc$antigenf)))
# identify the optimal change based on % agreement, take minimum change
optchgroc <- gridroc %>%
  group_by(antigenf) %>%
  mutate(maxagr=max(agreement,na.rm=TRUE),
         maxchg=ifelse( abs(agreement-maxagr)<0.001,1,0)) %>%
  filter(maxchg==1) %>%
  arrange(antigenf,abchange) %>%
  filter(row_number()==1) %>%
  select(antigenf,optchg=abchange,optagr=agreement,optkap=kappa,optcut=abcutoff)
# merge in the optimal change values
gridroc <- left_join(gridroc,optchgroc,by="antigenf")  %>%
  group_by(antigenf)

Figure to summarize results

The bottom panel plots the estimated seropositivity cutoff as a function of change in MFI levels used to identify presumed seroconverters. The mean plus 3 SDs of the distribution of antibody levels before presumed seroconversion was used to estimate the cutoff. For a given cutoff, the middle and top panels summarize classification agreement with the ROC-based seropositiity cutoff based on Cohen’s Kappa (middle) and proportion agreement (top).

Labeled dashed lines indicate the minimum magnitude of change in antibody levels that led to highest agreement, along with the corresponding seropositivity cutoff values, Cohen’s Kappa, and agreement.

plot_optchange(results=gridroc,classlab="ROC-based classification")

Sensitivity analysis over fold-change values, mixture model-based reference

Search over fold-change values to determine the cutoff with best agreement with mixture model-based classification

dmix <- dl %>%
  filter(!is.na(mixcut)) %>%
  mutate(antigenf=droplevels(antigenf))
gridmix <- foreach(ab=levels(dmix$antigenf),.combine=rbind) %:%
  foreach(icut=seq(0.3,2.5,by=0.01),.combine=rbind) %dopar% {
    # downsample to first measurements who increased by the increment icut
    # these are "presumed unexposed" under this cutoff in change icut
    # estimate mean and SD among the presumed unexposed to use as a 
    # seropositivity cutoff
    dmu <- dmix %>% 
      ungroup() %>%
      select(antigenf,time,logmfi,logmfidiff) %>%
      filter(antigenf==ab & time=="A" & logmfidiff>icut)
    
    # if there are <2 obs, cannot estimate mean + SD, so skip
    if(nrow(dmu)>=2) {
      
      dmu <- dmu %>%
        summarize(mu=mean(logmfi),sd=sd(logmfi),ucut=mu+3*sd)
    
      # classify all observations using the derived cutoff
      di <- dmix %>%  
      filter(antigenf==ab) %>%
      mutate(posi=factor(ifelse(logmfi>dmu$ucut,1,0),levels=c(0,1)),
             posmix=factor(posmix,levels=c(0,1)))
    
      tabi <- as.vector(table(di$posmix,di$posi))
      names(tabi) <- c("mix0cuti0","mix1cuti0","mix0cuti1","mix1cuti1")
      kappai <- psych::cohen.kappa(table(di$posmix,di$posi))$kappa
      tabi <- data.frame(t(tabi))
      tabi <- tabi %>%
        mutate(Nobs=mix0cuti0+mix1cuti0+mix0cuti1+mix1cuti1,
             agreement=(mix0cuti0+mix1cuti1)/Nobs,
             kappa=kappai,
             antigenf=ab,
             abchange=icut,
             abcutoff=dmu$ucut) 
      
    } else{
      tabi <- data.frame(mix0cuti0=NA,mix1cuti0=NA,mix0cuti1=NA,mix1cuti1=NA,Nobs=NA,agreement=NA,kappa=NA,antigenf=ab,abchange=icut,abcutoff=NA)
      
      }
    }  
gridmix <- gridmix %>% mutate(antigenf=factor(antigenf,levels=levels(dmix$antigenf)))
# identify the optimal change based on % agreement, take minimum change
optchgmix <- gridmix %>%
  group_by(antigenf) %>%
  mutate(maxagr=max(agreement,na.rm=TRUE),
         maxchg=ifelse( abs(agreement-maxagr)<0.001,1,0)) %>%
  filter(maxchg==1) %>%
  # mutate(agr99=ifelse(agreement>0.99 & !is.na(agreement),1,0)) %>%
  # filter(agr99==1) %>%
  arrange(antigenf,abchange) %>%
  filter(row_number()==1) %>%
  select(antigenf,optchg=abchange,optagr=agreement,optkap=kappa,optcut=abcutoff)
gridmix <- left_join(gridmix,optchgmix,by="antigenf") %>%
  group_by(antigenf)

Figure to summarize results

The bottom panel plots the estimated seropositivity cutoff as a function of change in MFI levels used to identify presumed seroconverters. The mean plus 3 SDs of the distribution of antibody levels before presumed seroconversion was used to estimate the cutoff. For a given cutoff, the middle and top panels summarize classification agreement with the mixture model-based seropositiity cutoff based on Cohen’s Kappa (middle) and proportion agreement (top).

Labeled dashed lines indicate the minimum magnitude of change in antibody levels that led to highest agreement, along with the corresponding seropositivity cutoff values, Cohen’s Kappa, and agreement.

plot_optchange(results=gridmix,classlab="Mixture model-based classification")

Session Info

sessionInfo()

R version 3.5.0 (2018-04-23)
Platform: x86_64-apple-darwin15.6.0 (64-bit)
Running under: macOS High Sierra 10.13.6

Matrix products: default
BLAS: /System/Library/Frameworks/Accelerate.framework/Versions/A/Frameworks/vecLib.framework/Versions/A/libBLAS.dylib
LAPACK: /Library/Frameworks/R.framework/Versions/3.5/Resources/lib/libRlapack.dylib

locale:
[1] en_US.UTF-8/en_US.UTF-8/en_US.UTF-8/C/en_US.UTF-8/en_US.UTF-8

attached base packages:
[1] parallel  stats     graphics  grDevices utils     datasets  methods   base     

other attached packages:
 [1] bindrcpp_0.2.2    doParallel_1.0.11 iterators_1.0.9   foreach_1.4.4    
 [5] gridExtra_2.3     psych_1.8.4       forcats_0.3.0     stringr_1.3.1    
 [9] dplyr_0.7.4       purrr_0.2.4       readr_1.1.1       tidyr_0.8.0      
[13] tibble_1.4.2      ggplot2_3.0.0     tidyverse_1.2.1   here_0.1         

loaded via a namespace (and not attached):
 [1] reshape2_1.4.3   haven_1.1.1      lattice_0.20-35  colorspace_1.3-2 yaml_2.1.19     
 [6] rlang_0.2.2      pillar_1.2.2     foreign_0.8-70   glue_1.2.0       withr_2.1.2     
[11] modelr_0.1.2     readxl_1.1.0     bindr_0.1.1      plyr_1.8.4       munsell_0.5.0   
[16] gtable_0.2.0     cellranger_1.1.0 rvest_0.3.2      codetools_0.2-15 knitr_1.20      
[21] broom_0.4.4      Rcpp_0.12.18     scales_1.0.0     backports_1.1.2  jsonlite_1.5    
[26] mnormt_1.5-5     hms_0.4.2        stringi_1.2.2    grid_3.5.0       rprojroot_1.3-2 
[31] cli_1.0.0        tools_3.5.0      magrittr_1.5     lazyeval_0.2.1   crayon_1.3.4    
[36] pkgconfig_2.0.2  xml2_1.2.0       lubridate_1.7.4  assertthat_0.2.0 httr_1.3.1      
[41] rstudioapi_0.7   R6_2.2.2         nlme_3.1-137     compiler_3.5.0